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In Palatini f(R) gravity, the parameters of the Schwarzschild - de Sitter solution as well as 
the whole interior solutions of compact objects are expected to change when compared to general 
relativity. We solve the Palatini field equations numerically in the case of the models f(R) = R+aR 2 
and f(R) = R - fj, 4 /R, and using the equation of state of Fermi gas. We show how the density 
profiles and the prediction for the maximum masses of white dwarfs (the Chandrasekhar limit) 
and neutron stars (the Tolman-Oppenheimer-Volkoff limit) are altered, and thereby conclude that 
observations on compact stars may be used to exclude alternative gravity models. 
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I. INTRODUCTION 

In general, general relativity is regarded as an effec- 
tive theory of gravity, its breakdown being practically 
inevitable at the Planck scale. Whether the breakdown 
could show at lower scales has been a topical question for 
long. Classical extensions to Einstein's theory have not 
only been studied in the context of inflation and early 
universe, but also to account for dark energy and dark 
matter [1]. Moreover, any suggested alternative gravity 
model should agree with the Solar System constraints 
and other experimental astrophysical data [2] . 

The f(R) gravity models have been among the simplest 
and most popular extensions of the Einstcin-Hilbert ac- 
tion. Originally suggested as models of inflation [3], they 
have gained interest as gravitational dark energy during 
past few years [4]. The action 

S = ^ J d 4 x^gf(R) + S m (1) 

generalizes the gravitational sector to involve non-linear 
interactions in the Ricci scalar: setting f(R) = R — 2A 
corresponds to the canonical Einstein- Hilbert action with 
a cosmological constant A. 

In general relativity, the metric g is usually taken as 
the only free gravitational field. The case, where both the 
metric and the connection T are dynamical variables, was 
also studied very early after the formulation of general 
relativity [5] . The fundamental theorem of Ricmannian 
geometry states that on any (pseudo-)Riemannian mani- 
fold there is a unique torsion-free connection for the given 
metric g that satisfies metric compatibility, and that con- 
nection is the Levi-Civita (Christoffel) connection. The 
compatibility condition Vg = for the only traditionally 
(physically) allowed covariant derivative equals to de- 
manding that test (physical) particles follow "the short- 
est paths on the manifold" [6] . More formally, geodesies 
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(defined by the connection) coincide with the extremals 
of the length function (defined by the metric) on the man- 
ifold. Choosing the covariant derivative with respect to 
the Levi-Civita connection, the act of parallel transport is 
intuitive: it preserves the inner products, which includes 
the norm and the sense of orthogonality. 

While this all might sound reasonable, the harsh fact 
is that there is no unique way to compare two tensors if 
they are elements of different tangent spaces. In fact, it 
can be said that we do not know if the geodesies truly co- 
incide with the shortest paths on our spacctimc manifold. 
What we do know very accurately, by experiment [7], is 
that all local effects of gravitation disappear, which im- 
poses the condition 

r^(P ) = o (2) 

for the connection coefficients in the local Lorentz frame 
at any spacetime point Pq. In all, there seems to be 
no immediately clear, physical reason to select the Levi- 
Civita connection from all the possible connections on 
the manifold. 

The formulation, where both the metric and connec- 
tion are left as free variables, is referred to as Palatini 
formalism, while insisting on keeping the connection as 
the Levi-Civita connection is called metric formalism. 
Starting from varying the action, however, the Einstein- 
Hilbcrt action (with or without a cosmological constant) 
yields the connection to be the Levi-Civita connection, 
if the connection is torsion-free. In the case of the gen- 
eral f(R) action (1) this is not the case: the metric and 
Palatini formalisms give very different gravity out of the 
same action. 

Perhaps somewhat surprisingly, many metric f(R) 
models for dark energy, such as f(R) = R — ti 4 /R with 
/1 4 ~ 10~ 104 m~ 4 , are shown to contradict the Solar Sys- 
tem measurements [8] . These models fail to produce the 
experimentally bound exterior solution outside a Newto- 
nian star. The corresponding Palatini models, however, 
reduce to the limit of general relativity at the Solar Sys- 
tem and are perfectly consistent with the experiments [9] . 
Nevertheless it seems that the dark energy Palatini f(R) 
models fail to produce the observed large scale structure 
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in the universe [10]. 

Moreover, models with extra terms such as f(R) = R+ 
aR 2 are not restricted by weak field effects. Such terms 
are expected to arise in quantizing gravity. In this work 
we study whether and how the structure of stars with 
strong gravitational fields is changed in the framework of 
f(R) models in comparison with general relativity. While 
the literature on the microphysics of compact stars is vast 
(that is on the effects that kick in gradually as the state 
of matter changes at high energies), the issue of changing 
effective gravity has been rather neglected. 



that is, the equation (8) is a fourth order differential 
equation for the metric tensor. This feature allows for 
several vacuum solutions (Tp„ = 0). 

Taking the covariant derivative of the left hand side 
of equations (3) and (6) yields Bianchi identities, which 
is shown in detail in [11]. The continuity equation takes 
the same form in both metric and Palatini formalisms, 

V M T^ - 0, (9) 

where V M is the covariant derivative with respect the 
Levi-Civita connection. 



II. F(R) GRAVITY 



III. THE PALATINI SOURCE EQUATIONS 



The field equations of Palatini gravity are attained by 
varying the action (1) with respect to the metric and the 
connection, 



F(R)R tw - ^f(R)g^ = 8nT^ 



gF{R)<r\ = 0, 



(3) 
(4) 



respectively, where F = df/dR (in general relativity, 
F = 1, Eq. (3) identifies as the Einstein equation and 
Eq. (4) as metric compatibility). The covariant deriva- 
tive V p = V p (r) and the Ricci tensor R^ v = i? M „(r) 
(and thus the Ricci scalar via R = g^ v R^(T)) are func- 
tions of the free connection. From the latter equation, 
however, we see immediately that the connection is the 
Levi-Civita connection for the conformally related metric 



!%, = { fr}[h] 
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= I $ v } [g\ + ^g pa {g^d,F + g^F - g^d a F). 

In the metric formalism the connection is taken to be the 
Levi-Civita connection of the metric g. The field equa- 
tion is obtained by varying the action (1) with respect to 
the metric only, 

FRfiv - l -fg^ - VpV„F + g^UF = SnGT^, (6) 

where the Ricci tensor and the covariant derivative are 
now defined directly by the Levi-Civita connection. 

Taking the traces of the field equations (3) and (6) 
highlights the difference between the Palatini and metric 
gravities. In the Palatini case the trace equation relates 
the Ricci scalar to the energy momentum tensor alge- 
braically, 



FR - 2/ = 8ttT, 



(7) 



which is also the case in general relativity. In the metric 
formalism the trace equation is a second order differential 
equation for the Ricci scalar, 



3DF + FR- 2/ = 8ttT, 



(8) 



Interiors of stars in equilibrium can be well described 
by static perfect fluid [14], 

T£ = diag(p(r),p(r),p(r),p(r)). (10) 

In this case, the metric can be parameterized 

ds 2 = -e A ^dt 2 + e B ^dr 2 + r 2 dtl 2 . (11) 

Supplemented by an equation of state, p = p(p), there are 
three functions to be solved, namely the metric functions 
A(r) and B(r), and the spatial structure of the star p — 
p{r). 

The continuity equation (9) reads 



A 1 

P' = -y(p + f>). 



(12) 



The rr and tt components of the field equation (3) take 
the forms 



A' = 



1 ( 1 - e B e B 



B' 



1 + 7 



1 / 1 - e B e B 



a 

—8nrp-\ 

F r 



(13) 



Hnrp f ) , (14) 



1 + 7 \ r F r 

respectively, as is shown in [9] . The variables 
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(16) 
(17) 



depend on the first and second derivatives of F. The 
trace equation (7) yields 



F" 



= N lP " + N 2 p 12 , 
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Nip', 



(18) 



(19) 
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where 
Ni 

N 9 



1_&F&R („_dp 
FdRdf \ dp 

1 \d 2 F fdR\ 2 dFd 2 R 
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(20) 



dp 
dp 



(21) 



Using the continuity equation (12), its derivative, and 
the 99 component of the field equation (3) , 



A" 



-A'- 
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~F 



-A'B' 



~-(A'-B') 



16tt 



2 e B T ee + e B (f - RF) - 2F" 



3F'- 



{B'-A'--)F' 



2 F 

we end up with the equations 

ap' 2 +bp' + c 

1 



(22) 
(23) 



P 



N 2 p'- 



2 V p + p 



2(1 + 7) 



N 2 p'- 



1-e 1 



d 
r 
e 

2F 



(FR-f) 



8ne B .7(4^37) \ e B 

-—p+^^—j + ^iFR-f) 

7(2-37)1 p + p 



8ne l 
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where 

a = rNi 



l--{p + p)N x 



b = 2[l- [p + p)N x ], 

l-e B 
r 
and 

1 



-(P + P) 



87rre 
F 



2F 



d= -r A' + 2^ + 1. 



(24) 

(25) 
(26) 

(FR-f)]27) 
(28) 



The equation (23) can be solved by 

a ^ 



P = 



-b±Vb 2 -4ac 
la ' 



a = 0. 



(29) 



In general relativity, a = 0. The "+" branch of Eq. (29) 
gives general relativity at the limit a —> 



lim 



-b + Vb 2 - 4ac 
2a 



while the "-" branch diverges. 

When supplemented by an equation of state p — p(p), 
the equations (23) and (14) accompanied by the equa- 
tion (24) form a closed set of differential equations for p 
and B. 




FIG. 1: The energy density profile of a static, spherical object 
with the central density 10 17 kg/m 3 and the equation of state 
of radiation, p — p/3, in two different gravity models: f(R) = 
R - ^/R with p 4 = 10 



! ! m 4 (red solid line) and general 



relativity with a cosmological constant A 
dashed line). 
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A. The exterior solution in Palatini f(R) 

As mentioned above, the trace equation (7) relates the 
Ricci scalar to the trace of the energy-momentum tensor 
algebraically in the Palatini formalism. For instance in 
the case of f(R) = R — p A /R the trace equation yields 



R = -A-kT ± V(4ttT) 2 + 3p 4 , 



(30) 



while for f(R) = R + aR 2 the trace equation is the same 
as in general relativity, R = —8irT. By the exterior met- 
ric of the star we imply the vacuum solution = 0, 
the modified Schwarzschild metric. Since, in that case, 
R and thus / and F are constants, the exterior metric is 
easily and uniquely solved from (13) and (14): 



2M r 2 , 
s = 1- — -yA e// 



2M r 2 , 
B = -hi I 1 Y A e// 



where M. the Schwarzschild mass, and 



A 
F 



are constants. For the model f(R) = R — p 4 /R, 



(31) 



(32) 



(33) 



(34) 



Although in general one cannot conclude the nature of 
the interior solution (p ^ 0) from the trace equation 
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alone, it can be done in the case of p — p/3, which is 
the equation of state of radiation and ultrarelativistic 
matter, and which also has T = 0. The case is partic- 
ularly interesting since as the matter in stars gets more 
and more relativistic, the equation of state gets closer to 
the ultrarelativistic one. The equations (13) and (14) 



A' 



1 



B' = 



e B r 



8tt 
~F 



8tt 

T' 



p 



•A 



A 



eff 



(35) 



(36) 



are the same as in general relativity except for the con- 
stant coefficient F and the effective cosmological con- 
stant (33). In the case of f(R) = R - p*/R, F = 4/3 
and A e yy is given by Eq. (34). Hence, no matter how 
small the parameter /i 4 is, there is a visible difference to 
GR (Fig. (1)). In the case f(R) = R + aR 2 , F = 1 and 
A e f f = A and there is no difference to GR. As shown in 
the following, this is actually opposite to what happens 
before the limit of the equation of state of radiation is 
reached: the /i 4 case implies practically no changes to the 
structure of the stars while the a model yields a notable 
change already at relatively small parameter values. 



IV. FERMI STARS IN F(R) MODELS 
A. Fermi gas equation of state 

By the Pauli exclusion principle, the pressure of a gas 
of fermions is non-zero even at zero temperature. After 
burning all the nuclear fuel that sustains thermal pres- 
sure, a star begins to collapse due to gravity. Ultimately 
the degeneracy pressure of electrons or neutrons (or, hy- 
pothetically, quarks) takes over, stabilizing the star to a 
white dwarf or a neutron star (or a quark star), respec- 
tively. If the pressure of the matter states the star goes 
through while collapsing denser and denser never prevails 
over gravity, the star is believed to form a black hole. 

In the completely degenerate case, the fermions occupy 
the lowest momenta states up to the Fermi momentum 



Pf 



f 3h 3 p 
\ 87rm 



1/3 



(37) 



where po is the rest frame energy density of the fermions 
of mass m. The pressure and the energy density of the 
spin-| Fermi gas are easily integrated analytically [12], 
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(38) 
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FIG. 2: The equation of state parameter of Fermi gas w (Eq. 
(46)) as a function of the rest energy density po- 



where E = 



- 1 1 , K = and 




Eq. (38) has the asymptotic forms 
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CO, 



(40) 



(41) 



where x = pp/m. When only the first terms in these ex- 
pansions are retained, we are left with the familiar poly- 
tropes of non-rclativistic {x — > 0) and ultra-relativistic 
(x —>■ oo) Fermi gas 



respectively, with 
and 

K 2 



ts 5/3 
P = K lPo j 



IS 4 /3 

P = K 2 Po , 



2/3 



TT J 20m(pM) 5 / 3 



1/3 



he 



TT J 8(/iM) 4 / 3 ' 



(42) 



(43) 



(44) 



(45) 



One should read p = 1 and M = m for neutron stars, 
while typically for white dwarfs p = 2 (number of 
nucleons per one electron) and M — m p (the average 
nucleon mass) [13]. 



87T 



PF 



dpp 2 E = K [sinh(t) - t] 



(39) 



The equation of state parameter of the pure Fermi gas 
(that is, p including only the degenerate fermions — in 
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FIG. 3: The lower limit of the inequality (47), p 2 = (l~3w)p 
(and the total energy density p) as a function of the rest 
energy density po in the case of a pure degenerate neutron 
gas. 



white dwarfs the plasma is composed of electrons and 
nuclei) , 



p lsinh(t) -8sinh(|i) +3t 
p 3 sinh(t) — t 



(46) 



approaches 1/3 at high (rest) energy densities as shown 
in Fig. (2). As reasoned in Section (III A), we would 
expect to see a change in the star profiles even at very 
low parameter p 4 values if the equation of state is close 
to 1/3. The question is, how close? It turns out that the 
Fermi equation of state parameter (46) is close enough to 
1/3 for uninterestingly high values of p 4 . One way to see 
this is to compare the p 2 and T(= (3w — l)p) terms in 
Eq. (30): at low energies w « 1/3 and —T » p 2 (for 
any natural value of p 2 ), but as we reach higher energies 
and w — > 1/3, the two terms start to approach each other. 
If 



pi 2 > (1 - 3w)p, 



(47) 



definitely there should appear a change in the star pro- 
file. In Figure (3) we have plotted this lower limit on p 2 
in the case of the Fermi equation of state. Comparing 
the plot p 2 (po) to the plot of p(po) we see that the fea- 
ture w — > 1/3 does push the limiting p 2 down by several 
magnitudes, but nevertheless not even close to any phys- 
ically interesting values. The value p 4 DE = 10~ 104 m -4 
that yields dark energy is simply so tiny that we cannot 
see any difference with respect to general relativity in the 
case of compact stars. 



B. The polytropic divergence 

The viability of Palatini f(R) gravity has been ques- 
tioned by Barausse et al [15] because the factor F" in 



the field equations diverge at p — ► (the surface of the 
star) for polytropes p = Kp$ with T > 3/2. The Fermi 
equation of state (38) approaches the form (42) — that is 
r = 5/3 — at the edge of the star, so the divergence might 
seem worrying. While it truly inflicts a qualitative differ- 
ence between general relativity and Palatini f(R) theo- 
ries, one should inspect the situation quantitatively. The 
diverging term is given by Eq. (18), and due to Eqs. (24) 
and (23), it is proportional to 



N 2 p' 2 cxC f{R) (T)(p+p) 



2 d 2 p 



dp 2 



n 3-2r 



(48) 



at the limit p — > 0, where C 'f(R) 
in [16], in the case of f(R) = R 
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aR 2 the factor 



As is noticed 



T 



oca, (T«l/a) (49) 
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suppresses the divergence to very low densities. The den- 
sity scale at which the divergence becomes of order unity 
is 



{K 2 /a) l /^- 2T \ 



(50) 



PS,a 

The maximum value — p$, a ~ 10 6 kg/m 3 — that will 
be encountered in this paper is the case of a — 10~ 16 m 2 
and K = K\ (Eq. (44)) for neutrons. 
In the case of f(R) = R — p 4 /R, 

1 



C f = 



P 2 yW- 



1 



32tt 2 |T| 3 ' 



x T 2 + 3-T) 2 

m»M 2 ) 



1 



(51) 



where T = 4nT/p 2 , and the density scale at which the 
divergence becomes of order unity is 

PS,, ~ (K 2 /pT 1/2T - (52) 

In the case of p 4 = 10~ 104 m~ 4 and K = K x (Eq. (44)) 
for neutrons, ps,^ ~ 10~ 7 kg/m 3 . 

Formally the limit pp —> can be taken, and p as well 
as p and po approach smoothly zero. Physically, however, 
the Fermi gas equation of state breaks down already at a 
higher density inside stars. White dwarfs have a gas enve- 
lope and neutron stars an outer crust, where the matter is 
not degenerate. The polytropic divergence does not prac- 
tically restrict the f(R) model parameters much: if we 
insist the validity of Eq. (42) until p s < 10 3 ^ (the den- 
sity of water), the parameters of models f(R) — R—p 4 /R 
and f(R) = R + aR 2 must take values p 4 < 10~ 71 and 
a < 1. These upper bounds are so large that they are 
already excluded by observation [17]. 

In this paper we have used the Fermi equation of state, 
Eqs. (38) and (39), to track profiles of stars down to as 
low densities as the efficiency (accuracy vs. time) of nu- 
merics has admitted. The divergences arise as approx- 
imated by Eqs. (50) and (52). These densities are so 
small in the models we have chosen to present in this pa- 
per that we entitle ourselves to cut the star profiles well 
before the divergences arise. 
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C. Results: density profiles and stability curves 

In general relativity, we are used to the so-called mass 
parameter m(r) = r (l — e _B< - r ^) /2, which is nothing 
but a reparameterization of the metric function B{r). 
In terms of the mass parameter, the field equation (14) 
takes the form 



,(53) 



which in the case of general relativity separates to the 
simple form m(r) = J Q dr4irr 2 p(r) + r 3 A/6. The in- 
terior solution is smoothly connected to the exterior 
Schwarzschild solution by identifying M = m(rf) — 
7jA e ///6, where rj is the radius of the object. When we 
speak of the mass of a star, we mean the Schwarzschild 
mass. 

The total energy of a static, spherical object has to be 
defined by the Schwarzschild mass M, since it completely 
governs the geometry outside the matter (except for the 
(effective) cosmological constant). The equilibria for a 
number of particles are the configurations that extrem- 
ize the Schwarzschild mass, and these configurations are 
indeed what we get by solving the (modified) Einstein 
field equations [18]. 

Although abstruse, the equation (23) for dp[p(r)]/dr 
(accepting only the "proper" solution, which yields gen- 
eral relativity as a limit) and the equation (53) for 
dm(r)/dr are first order differential equations with re- 
spect to r, and the solution is uniquely defined by the 
initial values m(0) = (equals to setting B(0) = 0) and 
p(0) = pc like in GR. Thus the central density pc may be 
safely used as the parameter labeling the different equi- 
librium configurations — and, moreover, to serve together 
with m(0) (or B(0)) as the initial values in solving the 
field equations. At the same time, for instance, the to- 
tal particle number or the radius of the object do not 
uniquely define the solution. 

In Figure (4) we present, as examples of our white 
dwarfs, the rest energy density profiles and the behavior 
of the metric function B{r) in three different cases: in 
Palatini f(R) = R + aR 2 with a = 1CT 16 m 2 and 1CT 17 
m 2 , and in general relativity, all starting with the same 
central density (pc = 10 15 kg/m 3 ) but ending with dif- 
ferent values of the Schwarzschild mass (2.21, 1.43 and 
1.36 Solar mass, respectively). The larger a is, the less 
centrally peaked the profiles seem to become: in this re- 
gion, gravity is weaker compared to general relativity. 

The dependence of the Schwarzschild mass M on the 
central density pc for different f(R) models is presented 
in Figures (7) and (9) for white dwarfs and neutron 
stars, respectively. In general relativity, the first max- 
ima of the plots are known as the Chandrasekhar limit 
and the Tolman-Oppenheimer-Volkoff (TOV) limit, the 
maximum white dwarf mass and the maximum neutron 
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FIG. 4: The rest energy density profiles po(r) of white dwarfs 
with the same central density (10 15 kg/m 3 ) in three dif- 
ferent models of gravity: f(R) = R + aR 2 (blue curve 



10 



~ 16 m 2 , green curve a 



10~ m ) and general rel- 
ativity (black curve). The metric function B(r) is presented 
for the same models (solid line is for the interior solution, 
dashed line for the continued exterior (Schwarzschild) solu- 
tion). The Schwarzschild masses of the stars are 2.21, 1.43 
and 1.36Ms„ n , respectively. 



star mass, respectively. In general, the stationary points 
dM(pc) / dpc — are called turning-points, and they may 
signal a change from stability to instability, or vice versa. 

In the case of static, spherically symmetric solutions, 
the stability analysis is covered by examining small ra- 
dial perturbations of the equilibrium configurations. The 
equations for the frequencies u> of oscillations of the radial 
displacements £(r, t) — £(r)e lu)t , the acoustical modes, 
are derived in the case of general relativity for example 
in [18]. It can be shown that the turning points sig- 
nal changes in the stability of different acoustical modes. 
The frequencies of the acoustical modes form an infinite 
discrete sequence 



(54) 



and depending on the sequence of the turning points sta- 
bility (oj 2 > 0) or instability (u> 2 < 0) of the configuration 
is implied. The logic, which mode changes at which turn- 
ing point, is straightforward: As we increase the central 
density, the radius R cither decreases or increases de- 
pending on the dominating mode in question. In the 
former case, a mode with an even number of nods dom- 
inates and changes stability at the turning point, while 
in the latter case the mode that changes stability must 
have odd nods [14]. Which even/odd mode is in order to 
change follows from the property (54). 

The key result is that the changes in stability can 
be traced directly from the M versus r/ plot since 
dM/drt — dM/dpc -dpc/drf. for instance, if all modes 
are stable at lower central densities, a counterclockwise 
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FIG. 5: The rest energy density profiles po(r) of neutron 
stars with the same central density (I0 18 kg/m 3 ) in three 
different models of gravity: f(R) — R + aR 2 (magenta curve 
a = I0 -19 m 2 , red curve a — 10 -20 m 2 ) and general rela- 
tivity (black curve). The metric function B{r) is presented 
for the same models (solid line is for the interior solution, 
dashed line for the continued exterior (Schwarzschild) solu- 
tion). The Schwarzschild masses of the stars are 0.239, 0.554 
and 0.625-Msun, respectively. 



FIG. 6: The rest energy density profiles po(r) of neutron 
stars with the same central density (I0 20 kg/m 3 ) in three 
different models of gravity: f(R) = R + aR 2 (magenta curve 
a = I0~ 19 m 2 , red curve a — I0~ 20 m 2 ) and general rela- 
tivity (black curve). The metric function B(r) is presented 
for the same models (solid line is for the interior solution, 
dashed line for the continued exterior (Schwarzschild) solu- 
tion). The Schwarzschild masses of the stars are 0.419, 0.494 
and 0.406Ms„ n , respectively. 



bend in the M{rj) plot at a critical point manifests an 
onset of instability with increasing pc (cf. the Chan- 
drasekhar and TOV limits). The behavior of dM/drj in 
different f(R) models is shown in Figures (8) and (10). 

From the white dwarf curves, Fig. (7) and (8), we 
immediately take notice of the disappearance of Chan- 
drasekhar limit at a > 10~ 16 m 2 . For lower values of a 
the instability is set on at the Chandrasekhar limit, but 
it seems to be followed by a region of stable solutions. 
However, the neutronization of the plasma is expected to 
begin relatively soon after the "Chandrasekhar critical 
density" — or even before it, depending on the compo- 
sition of the star — and to be an even more important 
factor in the collapse of white dwarfs than the instability 
due to gravity. 

Our stability curves for white dwarfs were plotted using 
the equation of state of degenerate electrons alone, that 
is, not taking the onset of neutronization into account. It 
is doubtful that the stability regions beyond the neutron- 
ization threshold pc ~ 10 15 kg/m 3 would act as barriers 
against collapsing. Nonetheless, since the Chandrasekhar 
limit appears to be missing completely for a > 10 -16 m 2 , 
the instability is to be triggered by the neutronization 
alone at a considerably higher M compared to general 
relativity. For a ~ 10~ 17 m 2 , the stability region be- 
gins slightly before the neutronization threshold, and it 
is an open question whether condensing taking place be- 
tween it and the Chandrasekhar limit is enough to set 
off nuclear reactions (and the explosion of the star). If 
it were not the case, the white dwarf needs to keep ab- 
sorbing matter from e.g. a binary partner until the nu- 



clear reactions start taking place at a higher value of the 
Schwarzschild mass. 

The stability curve M(pc) for neutron stars shown in 
Fig. (9) goes through a drastic change between a = 10~ 19 
m 2 and 10~ 20 m 2 . The Tolman-Oppenheimer-Volkov 
limit exists for a < 10~ 20 m 2 , but it is located at a 
lower Schwarzschild mass and at a lower central density. 
Moreover, the instability region subsequent to the TOV 
turning point seems to be followed by a stability region 
unlike in general relativity. For a > 10~ 19 m 2 , all solu- 
tions are stable and the TOV limit does not exist. 

Qualitatively, the behavior of the stability curve in the 
case of neutron stars is not as easy to understand as in 
the case of white dwarfs. From the Figures (5) and (6) 
we perceive that the density profiles of our neutron stars 
arc the less centrally peaked with increasing a. However, 
one should remember that the Schwarzschild mass has 
not much to do with the integrated energy density: it is 
defined by the value of the metric function B(r) at the 
edge of the star r = r/. The behavior of B(r) is also 
presented in Fig. (5) and (6). In the latter figure, the 
B{r) curve seems to be acting rather strangely in the case 
a = 10~ 19 m 2 : the mass parameter becomes temporarily 
negative inside the star. This is a general feature that 
starts to take place at the higher densities the smaller a 
is, but never in the case of general relativity. Before the 
edge of the star the mass parameter becomes positive, 
and the Schwarzschild mass is never negative. 

Naturally, at some high enough densities we expect 
to encounter the situation where 2m(r) jr > 1 — the 
black hole. In general relativity, a neutron star with mass 
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FIG. 7: The Schwarzschild mass to central density curves of 
white dwarfs in different models of gravity: f(R) = R + aR 2 
(blue curve a = 10~ 16 m 2 , green curve a — 10~ 17 m 2 , cyan 
curve a = 10 m , magenta curve a = 10 m , red curve 
a = 10 -20 m 2 and general relativity (black curve). 



exceeding the TOV limit is predicted to collapse due to 
the gravitational instability to these densities. In the a 
models, however, the collapse seems to be disrupted by 
a stable region like in the case of white dwarfs. Thus 
to form a black hole, the neutron star needs to absorb a 
considerable amount of mass from a partner star. This 
result, nonetheless, is questionable: the equation of state 
— the Fermi equation of state of degenerate neutrons 
alone, throughout the star — is expected to be far from 
reality. 



FIG. 9: The Schwarzschild mass to central density curves of 
neutron stars in three different models of gravity: f(R) = 
R + aR 2 (magenta curve a = 10~ 19 m 2 , red curve a — 10~ 20 
m 2 and general relativity (black curve). 



V. CONCLUSIONS 

In this paper we have studied the interior spacetimes of 
white dwarfs and neutron stars in Palatini f(R) gravity 
theories. A change in the structure of compact stars when 
compared to general relativity is expected both in the 
case of increasing and decreasing function F = df /dR: 
in the former case, this is due to the high densities (and 
thus large R) inside compact objects, and in the latter 
case due to the fact that the equation of state of the 
matter approaches the equation of state of radiation (and 
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FIG. 8: The Schwarzschild mass to radius of the star curves of 
white dwarfs in different models of gravity: f(R) = R + aR 2 
(blue curve a = 10 -16 m 2 , green curve a — 10 -17 m 2 , cyan 
curve a — 10 -18 m 2 and general relativity (black curve). 



FIG. 10: The Schwarzschild mass to radius of the star curves 
of neutron stars in three different models of gravity: f(R) = 
R + aR 2 (magenta curve a = 10 -19 m 2 , red curve a — 10~ 20 
m 2 and general relativity (black curve). 
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thus yields small R). 

To make a first sketch, we have used the equation of 
state of a single species spin-| Fermi gas from the center 
to the edge of the star, which is rather a rude simplifica- 
tion and could be taken to the next level to account more 
realistically for the structure of white dwarfs and neutron 
stars. After defining the equation of state, we are able 
to solve numerically the field equations of the particu- 
lar f(R) model and thus track the density/pressure and 
spacctimc curvature profiles of the stars. In all Palatini 
models, the interior solution is smoothly connected to the 
exterior Schwarzschild - de Sitter solution. Our results 
show how starting from a certain central density leads to 
different Schwarzschild masses in different f(R) models. 
In this paper we have mainly concentrated on the model 
f(R) = R + aR 2 and shown how 1) the Chandrasckhar 
limit disappears if a > 10~ 16 m 2 and 2) the higher the 
parameter a is, the lower Schwarzschild mass and central 
density the Tolman-Oppcnhcimcr-Volkoff limit is shifted 
to, and it seems to disappear for a > 1CP 19 m 2 . In the 
case of the model f(R) = R — ^ 4 /R we found out that 
there is no difference to the predictions of general relativ- 
ity even at extremely high densities, unless the parameter 



[i A » 1 m 4 , that is unnaturally large. 

The stability curves of white dwarfs and neutron stars 
can be easily produced in different f(R) models, and we 
will go through a selection of models and use improved 
equations of state in an upcoming article. The results 
can be used to restrict the parameter space of f(R) mod- 
els due to observations on binary systems and supernova 
explosions. While this paper was being completed, two 
papers, [19] and [20], on relativistic stars in metric f(R) 
theory also appeared on the arXiv. 
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